Mon. Not. R. Astron. Soc. 000, [HS] (2011) Printed 10 January 2013 (MN I^TeX style file v2.2) 



Period doubling bifurcation and high-order resonances in 
RR Lyrae hydrodynamical models 



o 



6 



in 



O 



X 



Z. Kollath*, L. Molnar, R. Szabo 

Konkoly Observatory, Budapest, 1121, Konkoly Thege ut 13-17, Hungary 



2011 January 26. 



ABSTRACT 

We investigated period doubling, a well-known phenomenon in dynamical systems, for 
the first time in RR Lyrae models. These studies provide theoretical background for 
the recent discovery of period doubling in some Blazhko RR Lyrae stars with the Ke- 
pler space telescope. Since period doubling was observed only in Blazhko-modulated 
stars so far, the phenomenon can help in the understanding of the modulation as well. 
Utilising the Florida-Budapest turbulent convective hydrodynamical code, we identi- 
fied the phenomenon in radiative and convective models as well. A period-doubling 
cascade was also followed up to an eight-period solution confirming that the destabil- 
isation of the limit cycle is indeed the underlying phenomenon. 

Floquet stability roots were calculated to investigate the possible causes and oc- 
currences of the phenomenon. A two-dimensional diagnostic diagram was constructed 
to display the various resonances between the fundamental mode and the different 
overtones. Combining the two tools, we confirmed that the period-doubling instability 
is caused by a 9:2 resonance between the 9th overtone and the fundamental mode. 
Destabilisation of the limit cycle by a resonance of a high-order mode is possible be- 
cause the overtone is a strange mode. The resonance is found to be sufficiently strong 
enough to shift the period of overtone with up to 10 percent. Our investigations sug- 
gest that a more complex interplay of radial (and presumably non-radial) modes could 
happen in RR Lyrae stars that might have connections with the Blazhko effect as well. 
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1 INTRODUCTION 

Period doubling, a phenomenon often observed in dynami- 
cal systems, was also found in stellar models and actual pul- 
sating variables in the last decades. Period doubling (PD) 
means that the observed quantity of the system alternates 
between a high and low amplitude cycle. Dynamical sys- 
tems as the simple Rossler oscillator are usually capable of 
period doubling bifurcation, and through a series of bifur- 
cations called the Feigenbaum cascade can evolve to chaotic 
behaviour. The very definition of RV Tauri variables was 
originally the altern ation of deep and shallow minima, a 
clear sign of PD (e.g. Prest on et al]|l963l l. The phenomenon 
was reproduced by iFokint (il99J) in radiative stellar models 
for RV Tauri stars as well as bv lSaitou. Takeuti fc Tanakal 
198911 in one-zone stellar models. Models of W Vir variables 
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cade towards chaos, and c haotic pulsations were indeed iden- 



tified in semiregular stars jBuchler et al 
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was reported by iKiss fc Szatmarvl \200i) in the Mira star 
R Cyg. Mo dels of classical puls a tors a lso showed promis- 
ing results. iMoskalik fc BuchleJ (| 19901 ) searched for half- 
integer resonances both in Cepheid and RR Lyrae mod- 
els and indeed found P D in the former case (see also 
iBuchler fc Moskaliklll990l ). A 2:3 resonance between the fun- 
damental mode and the first overtone was also identified as 
a root cause. But neither period doubling nor suitable res- 
onances were found in RR Lyrae stars between the funda- 
mental or first o vertone and a ny higher modes up to the 
fourth overtone. lAikawal (|200ll ) also reported PD in very 
long-period, radiative Cepheid models but did not identify 
any underlying resonance. 

The signs of possible PD are the half- integer frequencies 
(HIFs) in the Fourier spectrum, in the form of (2n + l)/2 /o 
with respect to the /o main periodicity, sometimes called 
as subharmonics. Similar peaks were found in some vari- 
able white dwarfs as well, suggesting PD in PG1351-I-489 
iGoupil et al.l 19881) and even sig ns of four-period behaviour 
in G191-16 i Vauclair et al.|[l98g ). But because exact values 
usually differ slightly from half-integer values, they could 
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Figure 1. Bifurcated solutions of a radiative model (6500 K, 0.57 Mq, 40 L0). Period doubling is shown in the upper left, four-period 
variation on the lower left panel. The right-hand panel shows the same models in the radial-velocity-radius phase space, the black solid 
line is the double-period where the blue dashed line is the four-period solution. Both models were calculated with 100 mass shells and 
identical parameters except the value of the artificial viscosity. 



be genuine non-radial modes that resonate with one o f the 
higher amplitude modes fe.o. lO'Donoghue et al.iri992l ). 

Period doubling in RR Lyrae variables was com- 
pletely unexpected. As not ed above, model calculations of 
iMoskalik fc Buchleij (jl990l ) did not yield positive results. 
No observational sign was ever reported either, beyond 
some residual scatter around pulsation maxim a and min- 
ima (iKoIenberg et al ] |2006l . lJurcsik et al.l [200^ 1. With the 
typical pulsation periods spanning between 0.2-1 days, only 
multi-site campaigns have had real chance to discover period 
doubling, by covering consecutive pulsation cycles. But even 
an extensive observing run on RR Lyr did not resolve any 
additional frequencies beyond the pulsational, Blazhko and 
combination peaks or any other hints of period doubling 
l|Kolenberg et all l200d 2011). More interestingly, PD had 
not been ide ntified unambiguousl y in the CoRoT datasets 
yet, although iPoretti et al.l ( 2010f ) mention a weak peak at 
1.5/0 for the modulated RRab CoRoT 101128793. 

The first signs of alternating pulsation amplitudes were 
found in RR Lyr itself through the contin uous and precise 
phot ometry of the Kep ler space telescope fKo lenberg et al] 
l20r0 b. GiUil and et al.l l2010l. Szabo ct al. (20l3) confirmed 
the case of RR Lyr and reported the discove ry of period 
doubl ing in two other RR Lyrae stars while iBenko et al.l 
l|2010l ) added four more Kepler suspects to the list where 
weak HIF peaks are present in the frequency spectra. In- 
terestingly enough, all seven stars where PD is present or 
presumed are Blazhko variables. Period doubling is domi- 
nantly present only at given Blazhko phases and even then 
the magnitude of the PD amplitude variations differ from 
one modulation cycle to another. All these facts considered, 
it is not that surprising that ground-based observations did 
not succeed. 

ISzabo et al.l (|2010l ) gave detailed description of the Ke- 
pler observations. It useful to recapitulate the basic obser- 
vational facts regarding the period doubling in RR Lyrae 
stars: in the case of RR Lyr variation in subsequent maxima 
reaches O.'^l but it amounts to a few hundredth of magni- 
tude in the other stars. It is interesting to note that in some 



phases of the PD, especially for V808 Cyg, the minima are 
changing barely, however at other phases they show as large 
deviations from c ycle to cycle as the maxima (see Fig. 6 in 
ISzabo et aLll2010l ). The long cadence data hints the possibil- 
ity of more complex features, e.g. four-period struct ures but 
short cadence data will be needed to confirm it. Sza bo et al.l 
(|2010h also reported that the magnitude of the PD variation 
and the amplitudes of most of the HIF peaks in the Fourier 
spectra change in concordance. The highest HIF peak is the 
3/2/0 in all three cases suggesting a 3 :2 resonance, but such 
mechanism was excluded by Moskalik fc Buchlerl (|l990l ) for 
RR Lyrae stars. The amplitudes of consecutive peaks do 
not decrease smoothly: there is a bump around the 9/2/o 
peak, pointing towards the possibility of some higher or- 
der resonance. It was also suggested with hydrodynamical 
calculations that period doubling occurs because of a high- 
order resonance between the fundamental mode and the 9th 
overtone. In this paper we investigate the numerical calcula- 
tions in detail to prove that the 9:2 resonance is indeed the 
underlying mechanism behind period doubling. 



2 PERIOD DOUBLING IN RR LYRAE 
MODELS 

The period doubling phenomenon was identified in RR 
Lyrae model calculations by the Florida-Budapest turbu- 
lent convective hydrodynamical code. Throughout the cal- 
culations we used convective models with OPAL95 opacities 
{Z — 0.0001 metal content) and 150 mass zones. The inner 
boundary temperature was set to 2.5 million K. Neither the 
stability of the static model nor the stability of the limit cy- 
cle solution depends significantly on the depth of the model. 
Increasing the inner boundary temperature up to 10 million 
K results only a few percent shift in the growth rates and sta- 
bility roots. However, we preferred the lower inner boundary 
temperature as it provides better resolution in the models 
at the same number of zones. Fo r det ails of the cod e we re - 
fer to iKollath fc Buchleij ^200)^ and iKollath et all l|2002f ). 
The relaxation method, an effective way of reaching limit 
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cycle solutio ns llStellingweri|[T974l 'l. and the Floquet stabil- 
ity analysis (|Floquetl 18831 ) are implemented and were used 
extensively in the calculations as well. 

Period doubling was first encountered in the models in 
the summer of 2009. Because of the complete lack of obser- 
vational or theoretical evidence at the time, even with the 
space-based results of CoRoT, only low priority was given to 
the results. We recognised the importance of the PD models 
when t he first batch of Kepler data arrived (jKolenberg et al] 
l2010bl ). The initial hydrodyn amical results were published 
along with the observations in ISzabo et al.1 (|2010l ). 

A bifurcated limit cycle is shown in the upper left 
panel of Figure [1] Alternating maxima are pronounced, with 
0.10 7?Q difference between cycles, while the effect is much 
smaller, only 0.002 Rq at minima in this case. Period dou- 
bling occurs not just in convective but purely radiative hy- 
drodynamical models as well. When calculating the radiative 
limit cases of convective models that show period doubling 
and using sufficient number of zones and low enough arti- 
ficial viscosity (100 zones and c, < 1.9 in our case), the 
limit cycle becomes unstable and bifurcates to a double- 
period solution. If Cq is decreased further, the model bifur- 
cates to a four-period solution (lower left panel in Figure [1] 
the period-doubling cascade will be discussed in Section [T2|l . 
The prominent appearance of PD in radiative models raises 
the question why it was not reporte d before. In our own sur- 
vey fo r double mode pulsation in lSzabo. Kollath fc Buchlerl 
ll2004) . we were interested only in the transient behaviour of 
the amplitude evolution of the modes and did not investigate 
the limit cycles. 

Because period doubling in Cepheids was connected to 
the half-integer resonances, we set out to investigate them 
in the RR Lyrae models as well. A large grid of linear 
models were created to identify the occurrences of vari- 
ous resonances. As lower-o rder overtones were excluded by 
iMoskahk fc Buchleij l|l990l ). modes up to the 14th overtone 
were all calculated. But handling of the large number of pos- 
sible resonances with the fundamental mode required addi- 
tional tools. 



2.1 Diagnostic diagram 

Models with a given Po/Pk resonance spread out to a 
2D surface in the three-dimensional Tef f-M-L (effective 
temperature-mass-luminosity) space. Analysing multiple 
near-resonance regions is not straightforward in such cases. 
We found however that these surfaces are not twisted in 
all directions; a useful two-dimensional projection of the 3D 
space can be constructed where the resonant models are 
confined to a narrow region. With this projection one looks 
"edge-on" at the surfaces. The relation we used is T^ff vs. 
150M — L where the mass and luminosity values are in so- 
lar units. The exact value of the mass coefficient is possibly 
different for each resonances but since none of them reduces 
to a ID line, we used the value 150 for convenience. Our 
calculations cover a more extended range of masses of lumi- 
nosities compared to the physical range of stars, in order to 
follow all resonances throughout the 150M — L parameter. 

The regions for different eigenmodes and resonance val- 
ues are plotted in Figure [2l Occurrences of resonances with 
normal modes are monotonically changing whereas reso- 
nance regions of strange modes are curved. This difference 




Figure 2. The diagnostic diagram. Resonance regions are plotted 
in a projection of the original Tf,ff — M — L space, the x axis 
being 150M — L. Overtone and resonance orders are labelled in 
ofc; form where k is the overtone order and the I lower index 
is the Pq/P). = 1/2 value. The strange modes 089 and 099 are 
quite distinguishable from the normal modes. Different colours 
represent different resonance orders. 



comes from the same behaviour that is noted at the modal 
diagram too (see Section [2. 2^ : normal modes follow an or- 
dered and well separated series with fairly similar period 
ratios that explains the similar slopes while strange modes 
deviate and have changing slopes. 

There are a number of other resonances with different 
overtones close to the 9:2 resonance of the 9th overtone, 
two of which are lower order than the strange mode. We 
concentrate on excluding these two as the higher-order (fc > 
10) normal modes are even less likely to have any effect on 
the fundamental one. The 7:2 resonance region of the 6th 
overtone and the 5:2 resonance region of the 4th overtone 
are well separated and distinctly different from the strange 
mode. Therefore model sequences covering a wide enough 
range of 150M — L values could discriminate between the 
modes and select the one that causes the bifurcation. To 
achieve this, the stability properties of the models had to be 
investigated. 



2.2 Strange modes 

As mentioned earlier, no suitable resonances were found be- 
tween the pulsational mode and some overtone up to the 
fourth order in previous studies. Higher order radial over- 
tones were usually disregarded because they normally are 
heavily damped and thus have no measurable effect on the 
pulsation. In special cases though high overtone modes de- 
viate from this simple picture, in the form of strange modes. 

Strange modes were first identified in highly non- 
adiaba t ic sys tems such as luminous helium stars where 
IWoodI (Il976l ) found large differences between a di abatic 
and non-adiabatic mode sequences. ICox et al.1 (Il980t ) 
also reported additional non-adiabatic modes in mod- 
els with high luminosity-mass ratios (Lq/Mq ^ 10'*). 
The term "strange mode" was introduced by ICox et all 
(|l980l ) simply to name the additional mode among 
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Figure 3. Modal diagram of an RR Lyrae model sequence. Pe- 
riods are scaled with the fundamental mode. Black solid lines 
denote the adiabatic modes whereas blue dashed lines denote the 
non-adiabatic modes. The 9/2 period ratio is indicated with the 
horizontal dashed line. Deviations are most pronounced around 
this period ratio, indicating the presence of strange modes. 



the lower pulsational modes. It was originally thought 
that strange modes require strong nonadiabaticity, but 
they were also identified later in weakly non-adiabati c 
Cepheid and RR Lyrae models (|Buchler fc KollathI [2OO1I ). 
Strange modes were ev entually found in classical pulsators 
ijBuchler. Yecko fc Kolla th 1997 ) and adiabatic models of 
massive stars ( Kiriakidis. Fricke fc Glatzell Il993l ) as well. 



Both studies concluded that strange modes are essentially 
acoustic waves that are trapped by a barrier in the form of a 
large sound speed gradient inside the star. Partial ionisation 
zones can act as very efficient boundaries that effectively 
decouple the star into two separate regions, both having 
their own oscillation spectra. Strange modes are therefore 
modes that are trapped between the barrier and the stel- 
lar surface and usually have peculiarly high growth-rates. 
Although strange modes can be identified in adiabatic mod- 
els as well, the differences are much more emphasized in 
the non-adiabatic models. Despite the extensive theoreti- 
cal work however, there is no unique and precise defini- 
tion yet. For more details on strange modes we refer to 
ISaio. Baker fc Gautsch\^ (j 19981 ). 

Strange modes can be identified on the modal diagram 
of our RR Lyrae model sequence in Figure [3] Adiabatic 
(black solid lines) and non-adiabatic (blue dashed lines) 
align well up to the sixth overtone. But around 9/2/o (or 
4.5 Po/Pk on the figure) modes appear without evident adi- 
abatic pairs, switching from overtone seven to ten at higher 
effective temperatures through avoided crossings. The ninth 
and tenth overtones can have large linear growth rates too 
(Figure|4|) but positive values are largely above the blue edge 
of the RRab instability region and mostly outside the iden- 
tified PD region. The overtone r egion 8-10 was foun d to be 
close to the instability limit by iGlasner fc BuchleJ (|l993l ) 
too. After the identificati on of strange mode s in classical 
pulsators, it was shown bv lBuchler et all (|l997l ) that the be- 
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Figure 4. Linear growth rates of radial modes. Starting with 
the second overtone, modes become more and more damped, but 
then overtone eight gets close to the unstable regime, indicating 
its strangeness. Overtone nine and ten could become excited, as 
well. 



haviour of growth rates can be explained by the occurrence 
of strange modes. We observe the same excursion towards 
instability too, as illustrated in Figure |4] But we are less in- 
terested in linear growth rates than the effects of nonlinear 
interaction between the strange mode and the fundamental 
mode. 



3 STABILITY OF THE LIMIT CYCLE 

The stability of the limit cycle and the presence of a destabil- 
ising strange mode is described by the Floquct analysis. For 
the required calculations we r efer to iStellingwcrf (1974) and 
iBuchler. Moskalik fc Kovacd |l99j). We give only a brief 
summary of the method here: after reaching a limit cycle so- 
lution, we perturb each independent variable in the model. 
Then the evolution is followed for one period and the dif- 
ference from the unperturbed limit cycle is approximated 
through a Jacobi matrix. The stability of the limit cycle is 
then given by the eigenvalues of the matrix in the form of 
Fk = exp(Aft + jfjife) where Fk are the feth order Floquet coef- 
ficients and Afc and are the corresponding exponents and 
phases. For the stability analysis of the limit cycle, we used 
our turbulent convective code. Then the variables in the Flo- 
quet matrix are the radius, velocity, turbulent energy and 
total energy of the zones. However, the results obtained by 
radiative models (where turbulent energy is ignored) do not 
differ significantly from the turbulent ones, as the additional 
eigenvalues related to turbulent energy do not play a role in 
the destabilization of the fundamental mode limit cycle. In 
general, the addition of turbulent energy in the linear sta- 
bility calculations introduces only damped thermal modes. 
Low-order vibrational modes can be associated with the cor- 
responding Floquet coefficients as (j)k ~ Po/A!(iiiod27r) and 
Afe ~ Kfc-Po, if the system is close to the blue edge of the 
instability strip. But mode identification becomes more and 
more complicated, especially for higher-order modes when 
the model gets farther from the blue edge, because <pk itself 
shifts and close to resonances the phase locks to tt (or zero) 
for any overtone. 
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iMoskalik k. BuchleJ l|l99d ) showed that in the vicinity 
of resonances the usually complex Floquet coefficients bifur- 
cate into a real pair, often creating a bubble in the corre- 
sponding exponents when plotted against Te//: the two real 
values first diverge then converge to the complex meeting 
point as the control parameter is varied. The corresponding 
phases discriminate between period ratios: for integer reso- 
nances {nPk « Po, n small and integer) the Floquet phase is 
(j)k = while for half-integer resonances {{2n + l) /2Pk ~ Po) 
the phase is (pk = n during the split of exponents. Period 
doubling occurs in the latter case when one of the exponents 
(the upper arm) becomes positive. 

Because of the more rapid variations of higher order co- 
efficients a dense grid was required in all control parameters. 
Model sequences with selected 150M — L values and with 
25-50 K steps in effective temperatures were calculated. We 
only considered the eddy viscosity (q:,^) from the convec- 
tive parameters as it does not change the equilibrium stellar 
structure but can have measurable effects on the stability 
of the limit cycle. At each temperature step was varied 
from 0.04 to 0.1 (or even from « to 0.1 in some cases). 
Such wide range in ai, is not expected in real stars but was 
required to identify the coefficients correctly. The step size 
was set to Aq^ = 2 ■ 10~*, which translates into 3-500 in- 
dividual models for each and every mass, luminosity and 
effective temperature value. Because the equilibrium stellar 
structure is not modified when ai, is changed slightly, the 
model can be relaxed from the previous limit cycle solution 
to the new one directly. Still, a large amount of CPU time 
was required to map the PD instability region. Our calcula- 
tions took several months on a 4-processor Intel Core2 Quad 
2.83 GIfz machine but a similar analysis using the turbulent 
convection strength parameter Qc would require far more as 
presumably every single model would have to be iterated 
close to the limit cycle first and then relajced in every step 
instead. 

3.1 Period doubling and resonances 

Figure [5] shows a simple map of the PD instability region 
plotted on the diagnostic diagram. The selection was based 
on the Floquet spectra: the plotted models exhibit a positive 
resonant exponent at = 0.05. At low 150A/ — L values 
the resonance regions of overtone six and nine overlap or are 
close to each other and the same is true for the fourth over- 
tone at higher values. The sixth overtone crosses the PD 
region but no positive exponents are observed at lower or 
higher temperatures and the PD region is almost completely 
separated from the fourth overtone. On the other hand, the 
non-monotone strip of the 9:2 resonance of the ninth over- 
tone clearly follows the PD region. It is interesting to note 
that the centre of the PD region is shifted to lower effective 
temperatures than the place of the actual linear resonance 
of the ninth overtone. We will discuss this shift in more de- 
tail in section 14.11 Still, we consider the 9:2 resonance the 
best explanation since the other resonances are much more 
separated from the margins of the PD region. 

3.2 Bifurcation cascade 

Successive period doubling events create a bifurcation cas- 
cade that could lead even to chaos. Our models show a pe- 
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Figure 5. Models where the resonant, positive Floquet exponent 
is observed at = 0.05, e.g. period doubling is present plotted 
over the diagnostic diagram (vertical blue columns). The models 
follow the 9:2 resonance of the 9th overtone. Only mean values of 
the resonance bands are plotted (grey dashed lines), except for the 
9th overtone where the edges are shown (thick black lines). The 
two lower-order normal modes and the the 8th overtone strange 
mode arc plotted with orange and black lines respectively. 

riod doubling cascade as well: Figure |6] displays the multi- 
plication of the number of different pulsation maxima ver- 
sus the Oi, eddy viscosity parameter. Towards the lower val- 
ues the model bifurcates to an eight-period solution through 
double- and four-period stages. The cascading we observed 
in both radiative and turbulent convective models confirms 
that a bifurcation cascade can arise in the limit cycle itself. 
Alternating amplitudes, resembling period doubling, might 
be g enerated through specific conditions as well (jPokinl 
[l99J)- This phenomenon however is a true bifurcation in 
a dynamical sense, which is confirmed by the cascade. We 
note again however that models with very low eddy viscos- 
ity parameters are interesting from the dynamical point of 
view: they do not represent actual stars. 

It is interesting to note in Figure[6]that period doubling 
occurs at two different ranges. The cascade is prominent 
on the left-hand side, at < 0.04. Then the model bifur- 
cates again around ~ 0.07. The behaviour of the limit 
cycle is in perfect agreement with the value of the corre- 
sponding Floquet exponent that drops from large positive 
values to negative with increasing a^, but grows back to 
small positive between 0.069 < < 0.097 thus generat- 
ing small-amplitude period doubling there. 



3.3 Examples of Floquet spectra 

The calculated Floquet spectra revealed surprisingly com- 
plex behaviour. Rapid variations in the exponents corre- 
sponding to higher order (fc > 5) modes required a fine 
step-rate in the ai, parameter as discussed above. Then the 
exponents were sorted: values from one edge were followed 
along the sequence to identify individual Afe(a^)) datasets. 
Two different selection criteria were set: any sets contain- 
ing <l)k ~ resonant phases were filtered out while sets 
containing ~ tt were differentiated from the remaining 
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Figure 6. Bifurcation cascade of an RR Lyrae model {T^ff = 
6500K,M = 0.59Mq,L = 56Lq). Each model correspond to a 
limit cycle solution with different eddy viscosity parameters. 
The plotted values are the maximum stellar radii of all different 
cycles. The upper panel shows the corresponding Floquet expo- 
nent. When the exponent is positive, the limit cycle is bifurcated. 
The vertical lines indicate the successive zero crossings of the ex- 
ponent. 



non-resonant exponents. With this method, no resonant ex- 
ponents were missed. If one calculates the models for a sin- 
gle value only, exponents related to resonances could be 
missed because the — n condition is valid only for a 
limited range of a,j in some cases. 

A relatively simple example is shown in Figure [T] The 
first top few dashed lines could be identified with the lowest- 
order modes. At lower values however, the exponents start 
to cross each other (they are separated in phase of course), 
making the further identification almost impossible: the 
ninth overtone for example is shifted to higher values than 
expected because of the resonance. The exponent we are 
interested in is plotted with black and orange lines: black 
means the phase-locked (real) solutions with (j>k = n (split 
into two arms) that are connected with intervals with com- 
plex Floquet exponents (orange lines). The value of the 
exponent varies rapidly, confirming that a single model 
with a badly chosen a^/ would indeed miss the unstable 
(Qi/ < 0.055) or even the real valued solutions. 

Not all Floquet spectra are that easy to explain: we 
encountered not just splitting into two arms but more com- 
plex structures suggesting some connection and/or interac- 
tion between individual resonances. This is especially true 
for models with higher effective temperatures where most 
of the exponents vary more rapidly. However, these findings 
do not affect our conclusions on period doubling and will be 
the scope of future research. 

We attempted to reconstruct the exponents that con- 
tain resonant, phase-locked portions in the vs. T^ff plane 
for a given oti, value, to identify various bubbles with the cor- 
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Figure 7. Floquet exponents versus the eddy viscosity parame- 
ter. Resonant values {(f)i^ = tt) are shown with thick black lines, 
the connecting resonant, not phase-locked intervals (t/n, < tt) are 
shown with orange solid lines. Dashed blue lines represent other 
exponents with < 0^ < tt. Exponents where 0j, = appeared 
were excluded for clarity. One coefficient becomes strongly unsta- 
ble for lower eddy viscosities (a,j < 0.055), creating PD in the 
limit cycle. Exponents at low values start to cross each other. 
Model parameters are: Te// = QiQOK, M = OSMq, L = SSL©. 

responding resonances. This turned out to be less straight- 
forward as we expected because of the various interactions 
between the exponents, making the identification along the 
T^ff values either difficult or ambiguous. We also note that 
the exact shapes of the bubbles depend strongly on the 
chosen ai, value as seen in Figure [7] Instead we created a 
schematic diagram (Figure [S]), showing the positions and 
approximate sizes of the bubbles, labelled with the most 
likely resonances. The diagram displays the strong presence 
of the two strange modes (089 and o9g) and the negative 
exponents of the resonances except of the upper arm of 
oQg, further confirming that the o9g resonance is the only 
possible explanation behind the period doubling. The high- 
order resonances with negative Floquet exponents do not 
have any significant effect on the pulsation, most proba- 
bly they cannot be observed. The highest order resonance 
known which results in observable distortion of the light- 
curve with a negative Floquet exponent is the 1:2 resonance 
of the fundamental mode and the 4th o vertone in s-Cepheids 
(jFeuchtinger. Buchler fc Kollathllioool ). 



4 DISCUSSION 

4.1 Shift from the linear resonance 

The identified PD region clearly shows some shift regarding 
the resonance region of the strange mode, as stated above. 
One has to bear in mind though that period ratios were 
calculated from linear eigenvalues: actual non-linear periods 
will differ from those values. The difference can be easily 
computed for the fundamental pulsational mode from the 
period of the limit cycle. Calculation of non-linear high- 
overtone periods cannot be done so easily, except for the 
PD region: because these are resonant models, the phase 
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Figure 8. Schematic diagram of the resonant Floquet exponents 
versus effective temperature. Model parameters are: 150M — L = 
22.5, ~ 0.06. The most likely overtones and resonances are 
indicated. The limit cycle is destabilised if the exponents turn 
positive. 

lock of the Floquet coefficients to tt guarantee that the pe- 
riod of the ninth overtone is exactly Po/4.5. In our set of 
models, the period of the fundamental mode is increased by 
« 0.002 d, approximately the same value for all parameter 
values. The non-linear period shift of the ninth overtone has 
the same order in the centre of the PD region but it has a 
higher variation with temperature. At the edge of the PD 
instability, the strange mode can shift by up to ~ 0.01 d or 
10% in respect of the linear value (Pg ~ 0.1 d) because of 
the resonant coupling, indicating the surprising strength of 
this particular 9:2 resonance. 

We have to note however that the properties of the 
strange mode or surface mode depends on the boundary con- 
ditions and zoning. Numerical uncertainties can also play a 
role in the shift of periods. The difference between equilib- 
rium and limit cycle models can also change the period of 
the strange mode. 

4.2 Effects of metallicity 

Metal content in the models was set to Z ^ 10"" but field 
RR Lyrae stars are sometimes more metal-abundant. We 
calculated linear models with Z = 0.001 and 0.004 values 
too. The appearance of the 9:2 resonance in the diagnostic 
diagram changed in two ways (see Figure [9l): first, it shifted 
to higher effective temperatures, secondly, it became more 
extended and curved. The gap in the Z — 0.004 band is 
not a plotting artefact but a real phenomenon. There are 
no resonant models for the intermediate 150M — L values 
(between about 40 and 50), even up T^ff ^ 72007^, the 
period ratio of the fundamental mode and the strange mode 
avoids the 9:2 value. 

One can readily assume that the PD instability region 
shifts to higher temperatures as well. If the PD region is 
about as narrow for higher metallicities as it is for Z = 
0.0001, observing PD in stars with independently measured 
effective temperature and metallicity values could provide a 
simple test of this mechanism. 
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Figure 9. Location of the 9:2 resonance of the ninth over- 
tone on the diagnostic diagram for different metallicities. Black: 
Z = 0.0001, blue: Z = 0.001, orange: Z = 0.004. Points are indi- 
vidual model calculations, curves represent the upper and lower 
envelopes of the T^ff (150M — L) projection of the surface of reso- 
nant models. Note that for Z = 0.004 the intermediate 150M — L 
values are devoid of the particular resonance. 

4.3 The case of RR Lyr 

Out of the three Kepler Blazhko-modulated stars showing 
the period doubling effect, only RR Lyr has effective temper- 
ature, luminosity and mass values published already. An ef- 
fective temperature of Teff = 6480± 120i(' was derived from 
multicolour photometry bv lSiegel (1982 ) . JHK near-infrared 
time-series photometry yielded the following parameters: 
Te/ / = 6330 ± 50K, L = 48 ± 5Lq and M = 0.63 ± O.OIMq , 
using a Z = 0.0012 metallicity value (ISoUima et al. 20081). 
By me ans of high-resolution spectroscopy. iKolenberg et al.l 
(|2010ah derived metallicity between Z = 0.001 - 0.003 and 
also calculated a luminosity of L = 49 ± 5Lq . The mass and 
luminosity values translate to a 150M — L* = 46 ± 6 value. 

RR Lyr would fit perfectly into our PD instability re- 
gion (Figure [Sjl if it had lower metallicity by an order of 
magnitude. Although we did not calculate nonlinear mod- 
els with higher metallicities, RR Lyr could also fit in the 
Z = 0.001 PD instability. Linear resonances occur between 
6700 — 6800K in this case. The calculated PD region has 
a width of ~ 300^^ below the linear resonance: extrapolat- 
ing this relation, RR Lyr would fall to the edge (^ 6400A") 
of the PD instability of higher metallicity. These estimates 
suggest the the higher effective temperature measurements 
{Teff > 6400A:) are preferred for RR Lyr from the PD point 
of view. We cannot however exclude numerical uncertainties 
as the cause of discrepancy between the position of RR Lyr 
and the PD region. 



5 CONCLUSIONS AND FUTURE WORK 

We presented the analysis of the period doubling phe- 
nomenon we encountered in RR Lyrae hydrodynamical mod- 
els. The behaviour of the models match the obs erved period 
doub ling in modulated Kepler RR Lyrae stars l|Szab6 et al.l 
[2O1O). Our investigations revealed the following: 
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(i) The underlying phenomenon is the destabihsation of 
the fundamental mode limit cycle. Not only period doubling 
was observed: a cascade up to eight-period solutions was also 
followed. 

(ii) The root cause behind the period doubling is a 9:2 res- 
onance between the fundamental mode and the ninth over- 
tone. This is the first example of such a high-order resonance 
that is able to destabilise the limit cycle. The ninth over- 
tone itself was found to be a strange mode. The resonance 
is strong enough to shift the period of the overtone as large 
as 10 percent. 

(iii) Floquet analysis was found to be a crucial (and 
also time-consuming) method to connect period doubling 
to the 9:2 resonance, together with our diagnostic diagram 
(150M - L vs. Teff). 

(iv) The strength of period doubling is sensitive to the 
convective parameters of the models. In this investigation 
we considered only the eddy viscosity parameter, there are 
however others, most notably the strength of the turbulent 
convective flux and /or the mixing length. 

(v) Increasing the metallicity shifts the resonance (and 
thus the PD instability) to considerably higher effective tem- 
peratures. If confirmed, this relation could help narrowing 
down the metal content of other RR Lyrae stars showing 
period doubling. 

There are about 30 known RR Lyrae stars in the Ke- 
pler field, half of which do not exhibit Blazhko modulation. 
No PD was found however in an y of the non-m odulated 
stars even below mmag level (Sza bo et all I2OI0I '). Yet the 
models have no problems producing period doubling in non- 
modulated stars. Different factors could explain this hiatus. 
For example the Kepler sample might be too small; HB evo- 
lutionary models would shed more light on the population 
of the PD instability region and on the time span the stars 
spend in and out of it; the modulated stars cover a wider 
parameter space that make the conditions required for PD 
easier to occur. 

The results suggest that PD can indeed play an impor- 
tant role in the understanding of the Blazhko effect. But 
is period doubling a cause or an effect? If the modulation 
is indeed caused by some kind of internal variation of th e 
stellar structure, like in the Stothers model (|Stothersll200^ . 
then PD might be another manifestation of this mechanism. 
Assuming that physica l parameters vary over the Blazhko 
cycle f see ISodon |2009| and references therein), these stars 
could sweep over a wider parameter range, providing better 
circumstances for the PD instability to occur. 

Since PD is a resonance phenomenon, resonances could 
also be connected to the origin of the Blazhko modulation 
(seeH ovac 3l2009l for a recent summary on the models) . Thus 
purely geometrical explanations as the magnetic oblique ro- 
tator model (Shibahashi 2000) are not well suited since there 
is no direct connection between the modulating process, the 
rotation, and the modal reson ances. On the other hand, the 
non-r adial resonance model (|Nowakowski fc Dziembowskil 
I2OOII ) links the modulation to a 1:1 resonance between the 
fundamental and a non-radial mode. As the PD instability 
demonstrates it, high-order modes can have a strong effect 
on the pulsation, suggesting that other resonances or some 
more complex interplay of radial and non-radial modes and 
resonances might occur in the stars. On the other hand. 



the same resonance that is responsible for the period dou- 
bling can result in different kinds of bifurcations. The am- 
plitude equation formalism provides possibilities to extend 
the hydrodynamical calculations to a more general study. 
These theoretical investigations may provide new theories 
for Blazhko phenomenon as stated in a companion paper 
(|Buchler fc Kollathll201ll ) . All these new developments point 
towards the possibility of a complicated mechanism behind 
the Blazhko effect that may include resonances or variations 
in the stellar structure or even the two processes simultane- 
ously. 

The fact that high-order resonances with strange modes 
play an important role in the observed behaviour of stellar 
pulsation arises lots of new questions. We have started ex- 
tended hydrodynamical surveys to test the effect of these 
resonances for other types of radial pulsators and for differ- 
ent parameter sets. The resonance of only two modes pro- 
vides period doubling, three-mode resonances however {e.g. 
a resonance among two low-order modes and a high-order 
strange mode) can also provide interesting effects like three- 
mode pulsation. Predictions for such modal interactions will 
be discussed in a forthcoming paper. 
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